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ERRATA 
July, 1978 

1. Pages i, l f 24, and 25 

The SCOPE operating system is SCOPE 2.1.3 not 3.4.4 

2. Page 4, Equation 1 should read: 

b + + o 

3x c* 8 V 8z c 2 


3. Page 7, line 18 should be corrected to read: 

however, a rectangular grid could not.... 

■ 

4* ’ Page 10 lines 3 and 4 should read: 

....consist of !5 rows and 4^ columns, while the lower one 
has £ columns and 4 rows. 

lines 12 and 13, " increasive" should read "increasing” 

5. Page 11, Figure 1.6 
Interchange "Row" and "Column". 

6. Page 12, line 7 should read: 

....points of the coarser.... 

line 18 should read: 

....which are not corner control points.... 

7. Page 18, equation 2 should read: ‘ , 

CU (Wy ° 7T) + W * + DU <ji u 

+ CL \\ • n) + TT • v £ + DL ^ « BET 

line 11 should read: 

...while w is the.... 

8. Page 19, line 13 should read: 

o = (w u - w £ ) • in 

9. Page 20, line 2 should read: 

....all panels lie within the Mach cone .... 

10. Page 31, lines 8, 10, and 12 
"1-1" should read "1=1" 




11. Page 32, line 6 should be corrected to read: 
point L are (ZM(1,L) * ZM(2,L) , ZM(3,L)} 

12. Page 35, line 13 should read 
.... If t!lROP^ "l , the User .... 

line 17 should reads 
. . . .Thus , if NACASE > 1 , . . . . 

13. Page 36, line 8 should be corrected to read: 

SUBROUTINE CCOF (CU, CL, TU, TL, DU, DL, NCT, NBIN ) 

line 11 should be corrected to reads 
SUBROUTINE CBET (BET, NBIN ) 

line 15 should be corrected to reads 
....BETA (81,.... 

delete lines 25 and 26s 

and the other quantities have been described previously. 

14. Page 39, common block in Figure 2.2 should be identical to 
common block on Fig. 2.10 

15. Page 43, line 24 should read: 

....that wake networks only.... 

16. Page 50, Comments It was not necessary to reset the 
composite boundary conditions for the wake networks. The 
edge matching procedure will override these boundary 
conditions for the wake networks and apply the correct 
conditions. 

17. Page 54, the description of NLOPT-9 should read: 

Difference between upper and lower surface perturbation 
potential.*. 

18. Page 55, the velocity parameters “v" are vector quantities 
"v" 

19. Page 58, line 1 should be corrected to read: 

..... . the side abut nothing 

line 3 should be corrected to read: 

....apply to network corner control points.) 

line 9 add the following sentence: 

The edge matching indicator is further explained in the 
program source code. 

20. Page 60 Isentropic pressure coefficient equation should 

read: Y 

C n = {[1 + 2~* M 2 (2u + u* + v 2 + w 2 )3^ - 1} (Isentropic) 

CO 
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21. Page 111 Equation on line 17 should read 

9 . s 

E «,.(<>, - ° 0 - x, - c v,) ! 

i=l - . J 


22. Page 113 under the heading "Panel", 3 should read 5 

23, Page 116 line 1 should reads 
....Along BF, these quadratics.... 

line 21 should reads 4 

.... edge being weighted 10 _ more.... 


User ' s Notest 


1. User's would be well advised to check the fidelity of the 
panel representation of the actual configuration. Unit 
normal vectors calculated from a flat panel representation 
do not necessarily converge to the correct value for the 
actual curved surface as the number of panels is 
increased. The program's interpretation can be checked by 
plotting MX vs X from the control point data illustrated in 
Figure 3.4. 

2. The discussion on page 12 on adjacent networks of differing 
paneling densities should be correct. The user is warned 
though that extensive check out of the technique has not been done. 


3. Recent experience in supersonic analysis indicates that the 
angle of attack and sideslip should be aligned with the 
compressibility axis for best results contrary to the 
discussion on page 31. 
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0.0 Summary 


The purpose of this report is to present sufficient instructions for 
running the advanced panel pilot code developed under contract 
NAS2-7729, "Development of a FLEXSTAB Computer Program" pertaining to 
Advanced Aerodynamics Technology (TASK III, Category 2 - Supersonic Steady 
Flow). This software was developed as a vehicle for numerical experimen- 
tation and it should not be construed to represent a finished "production" 
program. The code is based on a higher order panel method using linearly 
varying source and quadratically varying doublet distributions for compu- 
ting both linearized supersonic and subsonic flow over arbitrary wings and 
bodies. 


This User's Manual contains complete input and output descriptions. A brief 
description of the method is given as well as practical instructions for 
proper configurations modeling. Computed results are a^so included to demon- 
strate some of the capabilities of the pilot code. Results to date have 
shown very reasonable comparisons between the advanced panel pilot code and 
other theories and experimental data. Recent studies at The Boeing Company 
have indicated that configurations as complicated as a supercruiser fighter 
delivering a standoff weapon as shown in Figure 0.1 are within the capa- 
bilities of the pilot code. 

The computer program is written in FORTRAN IV for the SCOPE 3.4.4 operating 
system of the NASA-Ames CDC 7600 computer. The program uses overlay structure 
and thirteen disk files, and it requires approximately XSl 000 (Octal) cen- 
tral memory words. 





Figure 0 o l Supercruiser Delivering a Standoff Weapon 
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1.0 General Description 


The advanced panel pilot code is intended to solve a variety of boundary 
value problems in steady subsonic or supersonic inviscid flow. The 
typical situation is the analysis of the flow past a prescribed configur 
ation. The body surfaces are specified and the flow properties sought 
in a body- fixed (x,y,z) systen, called the g lobal coordinate system . The 
flow far upstream of the body is characterized by its Mach number M»and 
angles of attack and yaw, a and 6 respectively. 


Figure 1.1 Global Coordinate Syitem 

The resultant flow is considered to be a small perturbation of a uniform 
flow. The direction of this uniform flow relative to the (x,y,z) system 
is given by angles a. and B. . Compressibility axes (x r , y . z 1 are 
defined so that the x - axis is aligned with the undisturbed flow. Thus, 
the (x, y, z) system is produced by rotating the (x c> y c » z c ) system 
through B about the z axis, and then through a about the new y axis 


Figure 1.2 Compressibility Axis 


It is generally advisable to take the undisturbed flow direction to be 
that of the free-stream; i.e., to set cc = « c and 3 = B c . However, if 
the same configuration is to be analyzed at a fixed Mach number for a 
range of a and 3 , computing time can be considerably reduced by 

fixing a and B somewhere in the middle of the ranges of a and 
B of interest, and to regard the differences a - a c and 3 - B c as 
contributions to the perturbation flow. Results presented in Section 4.2 
show this can be satisfactory as long asot-a_and B-B_are less than one 

C t 

or two degrees for supe: onic solutions. Criteria for subsonic solutions 
have not yet been established but previous experience indicates differences 
as great as 10 degrees may be acceptable. 

If terms quadratic in the differences between the local velocity, pressure, 
and density and their values in the freestream are neglected, the flow is 
irrotational , with a perturbation velocity potential which satisfies, see 
ref. 1, 

jii + a&.o CD 

V 8 V 8z c 

where 


B 2 = 1 - M* 

CO 


The velocity is related to <t> by 


V = V + 

co 



in compressibility coordinates and 


V - v„ + ( 


3 $ 3 $ , 

3x * 3y * 3z ' 


in the global coordinates. 


On any solid surface, the basic boundary condition is that the total mass 
flux vector PV/Po, be parallel to the surface. To a first approximation, 
the perturbation mass flux vector is given by 


tr = ( B 2 51 |t ) 

w 3x 3y 1 


3z 


in compressibility coordinates. Thus the fundamental boundary condition 
is that 


(w + T m ) * n = 0 


on all solid surfaces (and any other stream surfaces), where n is the 
unit normal to the surface, directed into the flow. 



figure 1.3 Direction of Unit Normal Vector 

The method on which the advanced panel code is based is to represent the 
potential as a superposition of fundamental solutions of the partial 
differential equation (1). These fundamental solutions are sources and 
doublets whose locations are prescribed but whose strength must be 
determined so as to meet the boundary conditions. Specifically, they are 
distributed continuously (or nearly so) over all the wetted surfaces of 
the configuration under study, over whatever vortex sheets are shed from 
trailing edges, and, in special cases, over other surfaces as well, such 
as inlet faces, exhausts, etc. 
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The surfaces are approximated as a continuous network of quadrilateral 

d 

"panels" whose vertices are basic input to the program. The panels are 
not necessarily planar, but are divided into plane triangular subpanels. 

On each subpanel, the source strength - is approximated as linear in local 
planar coordinates, while the doublet strength is taken to be quadratic. 

The parameters of these linear and quadratic functions are put in terms 
of certain fundamental unknowns which include the source and doublet 
strengths at the center of each panel. These unknowns are then deter- 
mined by simultaneous solution of the algebraic system which results by 
stipulating that the boundary conditions be satisfied at a sufficient 
number of control points. Once this is accomplished, the potential and 
velocity fields are known. The pressure field can then be calculated 
from an appropriate pressure-velocity relationship, and forces and 
moments calculated by integration. 

1.1 PANELING 

In all but the simplest cases, it is useful to divide the configuration 
to be analyzed into a number of networks . These are the surfaces whose 
shape is the basic input to the program. For flow past a cone, a single 
network will suffice. It is simply the outer surface of the cone. 
Similarly, only a single network need be used to represent an isolated 
quadrilateral a wing of zero thickness. (However, in both cases additional 
networks are required to model wakes.) 

Usually the networks collectively describe the outer (wetted) surfaces of 
the configuration under study and its wakes. However, it is also possible 
to use networks which are not exposed to the fluid. For example, one 
way to model the flow past thin wings is to distribute sources on their 
upper and lower surfaces, and vortices or doublets over the mean (camber) 
surface within the wing. Also, to examine the flow at points off the body, 
one may put a network through those points without disturbing the flow 
simply by specifying zero jumps across the network in the potential and 
in the normal component of mass flux. 











1.1.1 Specification of a Single Network Surface . 

Every network surface is specified by giving the coordinates of an array 
of grid points which is basically quadrilateral. That is, the array 
consists of NM "rows" of grid points which each contain NN points, where 
NN is the number of "columns" of grid points. Often, as in the cases of 
flow past a cone or past a delta wing, a triangular array would be more 
convenient. This is easily accomplished by letting a single physical 
grid point belong to several rows or columns, as is illustrated in Fig. 1.4 
for the paneling of a triangular wing surface. Here we show three 
alternative arrangements of grid points, which are the intersections of 
the grid lines shown. In the case of paneling "A", if we decide (arbitr- 
arily) to label the lines parallel to BC as the "columns" of the network 
and the rays which meet at A the "rows", we have rows of 5 columns each, 
and point A belongs to all 7 rows. In paneling "B" point B belongs to 
5 columns and in "C", point C belongs to 7 rows. 

While the generally quadrilateral organization of grid points is suffi- 
ciently general for most purposes, some care must be taken. For example, 
however, a rectangular could not be used on a circular domain. However, 
it is easy to find a grid which is acceptable, as shown in Fig. 1.5. 

Naturally, not every arrangement of points is equally effective. If the 
triangles shown in Fig. 1.4 represent the surface of an arrow wing, for 
example, it would be desirable to locate the confluence of rows or columns 
at the forward vertex of the wing, since the conicity of the flow yields 
significant variations of the flow parameters in that neighborhood. See 
Section 4.3 for the results of using various panel ings to analyze a thin 
delta wing. 

The difference between a row and a column cannot be detected until the 
grid points are given numerical labels for purposes of identification. 

Grid points in the same column are always labelled consecutively, while 
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Col umns 


Col umns 


Figure 1.4 Alternative Paneling of an Arrow Wing 






the grid point labels of adjacent points in the same row differ by the 
number of columns. Thus, given the grid point labels shown in Fig. 1.6, 
the network on top would be said to consist of 4 rows and 5 columns, 
while the lower one has 4 columns and 5 rows. 

An important point to note is that the labeling of points within the 
network is used to distinguish the "upper" surface of the network from 
its "lower" surface. When the network is part of the surface of a solid 
body, the upper surface is often chosen to be the side exposed to the 
flow, while the lower surface is within the body. 

The sence of a network maybe defined by the vectors N and M. The vector 
N corresponds to a column of grid points directed in the direction of 
increasive grid point 1 abler, while the vector M corresponds to a row of 
grid points directed in the direction of increasive grid point labels. 

Then the vector N x M is directed out of the surface, or into the flow 
when the network bounds the flow as mentioned earlier. Thus, for 
either labeling or grid points shown in Fig. 1.6, the "upper" surface 
of the network is under the printed side of the page. 

1.1.2 Notes on Multi -Network Problems . 

Two points need to be kept in mind in dealing with more than one network. 
One is simply a matter of bookkeeping. All grid points on a configuration 
must be uniquely labeled; points on different networks must have different 
numbers. Further, grid points which are common to more than one network 
must have different labels when regarded, as grid points belonging to the 
different networks. 

The second point is that it is necessary to use exactly the same grid 
points for two adjacent networks with similar grid density on their 
juncture, unless one or both networks is free of doublets. The program 
will run if this requirement is ignored, but the doublet strength will 
then be discontinuous at the network boundary, which is usually disastrous 




in supersonic flow even if the resultant effective vortex system is weak. 

It is permissible for two adjacent doublet networks to have different 
grid densities so long as the grid points of one network are identical to 
a subset of the grid points of the other at their juncture. Moreover, 
the grid points of the finer network must lie precisely (within the 
accuracy of the machine for a supersonic solution) on the rectilinear 
curve (straight line segments) joining the grid points of the courses 
networks. Examples illustrating this point are shown in Fig. 1.7. Small 
mis-matches are allowed for subsonic solutions. 

A related point is that abutting doublet networks must abut along 
complete edges, i.e., their network corner gird points must concide. For 
example, at a wing-body juncture, the body surface must be split into 
six separate networks at axial lines and at cross-sections which pass 
through the points at which the wing's leading and trailing edges meet 
the body. That is, of the arrangement of networks shown in Fig. 1.8, 
the upper one is correct but the lower is not. In the lower one, the 
two wing networks (its upper and lower surfaces) meet the body at points 
which are not control points for the body surface panels; therefore, 
continuity of doublet strength cannot be enforced at those points. In 
the upper arrangement, the wing networks meet the body at points which 
are corner points of four body-surface networks, and which therefore are 
control points. 

1.1.3 Surface Fitting . 

The four grid points on two adjacent rows and two adjacent columns deter- 
mine a quadrilateral, which need not be planar. 
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Acceptable for subsonic solutions only 


Figure 1.7 Adjacent Networks with Different 
Paneling Densities 
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Figure 1.8 Splitting of Body Surface into Networks 
at Wing-bodv Juncture 
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Figure 1.9 Subpanei Geometry 

However, the midpoints of the sides of this quadrilateral (shown as 0's 
in Fig. 1.9) determine another quadrilateral, which is planar (in fact, 
it is a parallelogram). Therefore, each set of four adjacent grid points 
(i.e., points like the dots in Fig. 1.9 which lie on two adjacent rows 
and two adjacent columns) is used to generate an approximating surface 
which consists of four triangles and one parallelogram, the corners of 
the parallelogram being the midpoints of the lines joining the four grid 
points. For computational convenience, the parallelogram is further 
subdivided into four triangles as shown by the dotted lines in Fig. 1.9. 
Since they all lie in the same plane, this subdivision is without geo- 
metrical significance. However, as is described in Appendix A, different 
doublet distributions are used on the four triangles within the parallo- 
gram. 

The nonplanar surface generated by the four adjacent grid points is 
termed a panel . As described above, each panel consists of eight tri- 
angular (and hence planar) subpanels , four of which comprise a parallel- 
ogram. 

1.2 SINGULARITY DISTRIBUTION 

Generally, every panel of every network is covered with a continuous 
distribution of sources and doublets. The source distribution is linear 
in the two coordinates which describe position on the panel, while the 
doublet distribution is quadratic. As described in Appendix A, the 
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parameters which describe the variations of both the source and doublet 
strengths on the panels are related to certain basic unknowns, namely 
o the source strength at the center* of every panel; 
o the doublet strength at the center of every panel; 
o the doublet strength at the midpoint of every panel edge 
which lies on the boundary of the network; 

and 

o the doublet strength at every corner of the network. 

This is accomplished in such a way that the doublet strength is continuous 
between subpanels and between panels. However, neither the source 
strength nor the derivative of the doublet strength normal to the edge 
is exactly continuous between panels. 

Naturally, the number of equations used to determine the singularity 
strength must match the number of unknowns. These are generated by 
imposing appropriate boundary conditions at selected "control points". 
Specifically, two boundary conditions are imposed at every panel center, 
one at the midpoint of every panel edge which lies on the network 
boundary, and one at every network corner point**. Fig. 1-10. 

The number of control points on a quadrilateral network (as shown in 
the upper part of Fig. 1.10) is (NM + 1) x (NN + 1), where NM is the 
number of rows and NN the number of columns of the network. When one 


* By "center" of a panel, we mean just the average position 
of its four corner points. If the panel has only three 
distinct corner points, the grid point which is input twice 
is counted twice as heavily as the other two. 

**To avoid singularities in the perturbation velocities due to 
discontinuities in the value or slope of the singularity distri- 
butions at panel or subpanel boundaries, control points are 
actually displaced slightly from the geometric center of the panel 
(which is the vertex of subpanels) and from the network edge. 
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edge of the network collapses to a point as in the lower part of 
Fig. 1.10, this number is reduced automatically. 

What we are describing here is the usual case in which the network is 
to be covered with both sources and doublets. The user has the option 
of suppressing either singularity type, as will be described immediately. 
Other options are available for the treatment of wakes. 

1.3 BOUNDARY CONDITIONS 

The general form of the boundary condition imposed at a control point 
is: 

CU (w u • "n) + TU 0 v u + DU <j> n 

+ CL (vJ^ •'n) + TT * + DL <{>£ « BET 

Here CU, TU, DU, CL, TL, DL and BET are to be specified, while i’s the 

perturbation mass flux vector, v the perturbation velocity vector, and <f> 
the perturbation potential. The subscripts u and 1 distinguish upper 
and lower surfaces of the network, respectively, and n is the unit normal 
directed from the lower to the upper surface of the network (i.e., in 
the NxM direction) . TU and TL are tangetial vectors. While the code 
contains both design and analysis type networks, only the use of analysis 
networks is discussed herein. 

At every panel-center control point, two boundary conditions of Eq. (2) 
must be specified, unless the network is specified as source- or 
doublet-free in which case only one condition like Eq. (2) is required. 
For control points on the edge or corner of a network, only one boundary 
condition like Eq. (2) is imposed, unless the network is doublet-free, 
in which case none are required. 


1.3.1 Boundary Conditions on Closed Impermeable Surfaces 
Under rather loose conditions, the potential at any point in the flow 
may be represented by distributing sources and doublets over the surfaces 
which bound the flow. The contributions of the sources and doublets 
may be split in various ways more or less according to the taste of the 
analyst. Some analyses are based almost entirely on source distributions, 
while others rely solely on doublets. However, Morino et al (ref. 2) 
showed that numerical problems in determining the singularity strengths 
were reduced if they were, in effect, rendered unique by stipulating 
that the perturbation potential vanish inside the body. So doing has 
the following effects: 

1. The source strength 

a m GT U ~ w £ ) • n 

becomes identical to the normal component of the perturbation 
mass flux on the outer surface of the body, w y - n which, if 
the surface is impenetrable, is simply the negative of the 
normal component of the freestream velocity. In other words, 
the source strength is known from the condition of impenetra- 
bility to be 

cr = - V* 77 

CO 

2. The doublet strength reduces to the perturbation potential on 
the outer surface, <f> u . Hence the perturbation velocity on 
the body surface can be computed by simply differentiating the 
doublet strength. 

In the pilot code, Morino-type boundary conditions are imposed by setting 
Cw u “ w^) *77= -V„.n at all panel-center control points, and ^ = 0 

at both panel-center and network-edge control points. 
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In supersonic flow these boundary conditions may only be applied on 
subinclined networks (i.e., all panels lie with the Mach cone from 
their respective leading edges). Superinclined boundary conditions are 
discussed later in Section 1.3.3. 

1.3.2 Hake Boundary Conditions. 

All doublet or composite networks must be followed downstream by another 
doublet or composite network or by a wake type network. This is necessary 
even for supersonic trailing edges to insure correct edge matching. A 
pure doublet wake network may be one of two types (designated "type 18" 
or "type 20"). Fig. 1-11 illustrates the use of these two network types. 

A wake network of type 18 is used behind a wing, with the boundary 
condition w- TT = - • rT specified on the leading edge, and no other 

boundary conditions being imposed. (This network type has only leading 
edge control points; there are no other control points.) The doublet 
strength on a type 18 network varies in the spanwise direction, but is 
constant in the streamwise direction. This enables the kutta condition to 
be satisfied. A type 18 network is also used behind the body in insure 
that the flow leaves the surface tangetially. 

A wake network of type 20 is used behind a curved body, and is a surface 
of constant doublet strength. The boundary conditions imposed are the 
same as before v7 • n = - V 03 *n’ but this condition is applied at the origin 
of the M-N coordinate system only, on the leading edge control point. 

Again there are no other control points. Often the leading edge of such 
a network will be collapsed to a point, this point being the intersection 
of the wing trailing edge and the body. In both types of wake networks, 
the distinction between rows and columns is no longer arbitrary. The 
M vector must always point in the down stream direction as shown in Fig. 
1.11. A type 20 network must not attach along a supersonic edge-o.k. , if 
wake network attaches at a point. 







1.3.3 Superincl ined Boundary Conditions - Supersonic Flow Only. 

In the analysis of the complete aircraft configurations in supersonic 
flow, one has to take account of the flow into the inlets of the jet 
engines. Where primary interest is in the pressure distribution on the 
exterior surfaces, the nacelle may be effectively closed by using a 
network of superincl ined panels (panels inclined ahead of the Mach 
cone from their respective leading edges) to cover the inlet. Boundary 
conditions can be applied only to the downstream side of any superincl ined 
network. In subroutine INPUT a superinclined network is designated as 
a composite network in the same way as for subinclined networks, see 
Fig. 1.10. With the downstream side denoted by the subscript 1, the 
boundary conditions 


<f> 


5L 


= 0 


are applied, the former at panel center control points and the later 
at both panel -center and network edge control points. These boundary 
conditions produce free stream flow downstream of the network. Again, 
note that superinclined boundary conditions are only applied in super- 
sonic flow. 


1.4 FORCES AND MOMENTS 

Forces and moments on each network are calculated in the following manner. 
Let S be a portion of a physical surface. The force coefficient vector Cp 
on S is defined as 

W // c P " ds 

where S is a specified condi guration reference area and h is the unit 
normal to the surface S (in the direction N x M). 
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A pressure coefficient can be defined for each side of S since there is 
always flow on both sides of S (although the flow on one side may be of 
no physical interest, e.g. inside a thick wing). The side of S on which 
the normal n points outward is called the upper surface and the other 
side of S is called the lower surface. Hence, both an upper and lower 
surface force coefficient can be defined. If the flow on both sides of 
S is physically significant (e.g., on a thin wing) then the total force 
coefficient for S must be considered as the sum of the upper and lower 
surface force coefficients. An upper surface, lower surface, and “sum" 
force coefficient are computed for each network. Moreover, these coeffi 
cients are accumulated separately over all (non-wake) type networks. 

The accumulated values are useful only when the physically significant 
force coefficient happens to be upper surface, lower surface, or "sum" 
for all networks of a given configuration. For example, the total 
force coefficient for a body and thin wing combination would have to be 
computed by hand by adding the upper surface body force coefficient to 
the "sum" wing force coefficient. 

The moment coefficient at a point about an axis t^ is defined by 

= t R * II S ^ ’ "V x n C p ds 

where T^ is a configuration reference length for the axis t^. Moment 
coefficients are calculated for t^ equal to each of the three coordinate 
axis. An upper surface, lower surface, and "sum" value is computed for 
each of these coefficients and these are to be interpreted in the same 
manner as described for the force coefficients. 



1.5 ; PROGRAM DESCRIPTION 

The computer program is written in FORTRAN IV for the SCOPE 3.4.4 operating 
system of the Ames CDC 7600 computer. The program uses overlay structure 
and thirteen disk files (including input and output files), and it requires 
approximately 132000 (Octal) central memory words. This program has been 
designed with the primary objective of verifying new concepts and numeri- 
cal experimentation. 

The computer program consists of one main overlay, eleven primary overlays, 
and two use libraries. The program flowchart shown in Figure 1.12 illu- 
strates the overall functions the program performs. 
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2.0 Input Data Description 

All input quantities are fed into the program in the subroutine called 
INPUT. A variety of INPUT routines were used during the course of the 
Category 2 work. It was found that in the majority of cases it was more 
efficient to use FORTRAN coded input than data card input. This is 
because about 90% of the test configurations were regular in shape and 
corner point data could be easily generated with code. On those 
occasions when externally generated corner point data v/as employed simple 
read statements were inserted at appropriate places in the code. 

As is summarized in Fig. 2-1: the input sequence is divided into four 
basic sections: 

1. The General Specifications section defines the basic parameters that 
control the amount of printout, defines the number of solutions, 
defines the reference parameters, etc. 

2. The Flow Conditions seciton defines the combinations of angle of 
attack and sideslip to be analyzed. 

3. The Network Geometry Specificaitons define the network geometry. This 
can be generated in the INPUT subroutine by the appropriate FORTRAN 
statements or generated externally to the program and read in on cards. 

4. The Boundary Condition Specifications section tends to be the most 
complex. Boundary conditions for each panel center and panel edge 
control point are specified in this section as well as the specifica- 
tion of network edge matching conditions. 

These sections will now be described in detail. An alphabetized glossary 
of the input variables is given in Table 2.1. 

2.1 GENERAL SPECIFICATIONS 

Variables related to the basic problem specification which must be 
supplied in INPUT include: 
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GENERAL 

SPECIFICATIONS 


Set: IGE0MP, IC0NTP, ISINGP, ISINGS. IBC0NP , IEDGEP, IPRAIC, 
AMACH ,ALPC , BETC , XREF , YREF, ZREF, SREF. BREF, CREF, 
DREF , NS YHH , NN ETT , NACASE .NDTCHK 


FLOW 

CONDITIONS 


For 1=1 to NACASE 
Set: ALPHA(I),BETA(I),FSVM(I) 


Set: L“0,ICA=0 (cumulative grid and* panel indices) 


NETWORK 

GEOMETRY 

SPECIFICATIONS 


BOUNOARY 

CONDITION 

SPECIFICATIONS 


For K=1 to NNETT 

Set: NTSfKLNTDfKLlP^TpO ,NM(K),NN(KK) 



For IC=1 to (NN(K)+1}*(NM(K)+1) 
or (NN(K)+1) 
or 1 

Set: ICA=ICA+1 

Set: NL0PT1 

NR0PT1 (IF HL0PT1^O) 
CU!,CL1,TU1,TL1,DU1,DL1,NCT1(IF NL0PTl=l) 
BET(I) ,1=1, NACASE (IF NR0PTl=l) 

NL0PT2 

NR0PT2 (1 F NLfiPT2i E 0) 

CU2.CL2 ,TU2 ,TL2,DU2,DL2 ,NCT2( IF NL0PT2=1 ) 
BET(I) , I=1,NACASE(IF NR0PT?=1) 

Call: IBTRNStICA.CUlKto store/BCON/on tape 1) 


Figure 2.1 Input Data Sequence 
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- 4" V: K : <r:rr.. 


AMACH 

ALPC, BETC 
NACASE 


NSYMM 


NNETT 

SREF 


= Mach number of flow far upstream (M© 

= Angles defining orientation of compressibility 
axes a c , $ c relative to body- fixed axes. 

= Number of different combinations of a and jS 
to be studied for the same configuration, a c s 
3 c , and M^, (must not exceed 8). 

= 0 if the flow is not symmetrical about the 
x-z plane. 

= 1 if the flow is symmetrical about the x-z plane. 

= 2 if the flow is symmetrical about the x-z and 
x-y plane. 

= Number of networks (must not exceed 50). 

- Reference area for dimensionless force and 
moment coefficients. 


(XREF, YREF, ZREF) = Components of point about which moment coefficients 

are computed. 

(BREF, CREF, DREF) = Reference lengths for computation of (x,y,z) 

components of dimensionless moment coefficient. 


The pilot code presently has the following size limitations: 
1 4 NNETT 50 
1 « NACASE « 8 


In picking NNETT, the necessity for wake networks should not be forgotten; 
see Section 1.3.2. 
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There are default options for most of the above-listed parameters; see 
Table 2.1. In particular, note that, unless otherwise stated, the 
reference lengths for the moment computation are each 1.0, and the moments 
are taken about e origin. These default options are often satisfactory. 
The amount of output may be controlled by variables set in INPUT. See 
Section 3 for a description of the output. 

IGEOMP = 1 if geometry diagnostic data (in addition to 

mesh point data) are desired* 

ISINGP = 1 if singularity spline data are desired’* 

ICONTP — 1 if control point and edge downv/ash data are 

desired 


IBCONP = 1 of boundary condition defining parameters are 

desired 

= 1 if edge matching diagnostic data are desired 

= 1 if singularity strength data are desired 

= L if influence coefficient diagnostic data is 
desired for the Lth control point 


IEDGEP 

ISINGS 

IPRAIC 


NDTCHK = 1 for program stop immediately before commencement 

of influence coefficient calculations. This can be 
used as a data check before committing the program 
to full execution. 


The default value of these variables is 1, so that the corresponding 
data are output unless they are set equal to 0 (except for IPRAIC, whose 
default value is 0). 


*Usually these data are not required. 
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2.2 FLOW CONDITIONS 

As many as eight combinations of angle of attack, sideslip, and free- 
stream velocity magnitude may be specified in a single run. Results may 
be suspect if the specified angles of attack and sideslip differ from the 
compressibility axis angles by more than one or two degrees in supersonic 
flow. 

ALPHA (I) 

BETA (I) 

fsvm (i) 

2.3 NETWORK GEOMETRY SPECIFICATIONS 

The panel corner point geometry may be generated within the INPUT subroutine 
by appropriate FORTRAN code or the geometry may be generated by some 
external paneling programs and read in on cards. Sometimes, especially 
when the network geometry is given as the punched-card output of some pre- 
processor program, the grid points are ordered in such a way as to make 
the outward normal point in the wrong direction. When this is the case, 
a subroutine MNSWCH may be called to make column into rows and vice versa: 

CALL MNSWCH (ZM( I, L) , NM(K) , NN(K)) 

Here K is the identifier of the network whose rows and columns are to be 
switched, and L is the label of the first grid point in that network. 

MNSWCH also interchanges the row and column numbers NM(K) and NN(K) in 
the process. 


= Angle of attack at upstream infinity relative 
to body-fixed axes (I-l, , NACASE) 

= Sideslip angle at upstream infinity relative 
to body- fixed axes (I-l, NACASE) 

= Freestream velocity magnitude-Default value 1.0 
(I-l NACASE) 
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The input variables defining the networks are: 


NM(K) = Number of rows in Kth network (must be ^ 2) 

MN(K) = Number of columns in Kth network (must be ^ 2) 

ZM(I S L) = Ith coordinate of Lth grid point (1=1, 2, 3) 

That is, the global coordinates of the grid 
point L are (ZM(L,L), ZM(2,L), AM(3,L)} 

Also, the pilot code is currently limited to a total of 1500 grid points, 
1500 control points and 1000 panels, which may be distributed at will 
among the different networks. 

In setting up neighboring networks, the points raised in Section 1.1.2 
must be kept in mind. 

2.4 BOUNDARY CONDITIONS 

One or two boundary conditions must be imposed at every control point, 
depending on the type of network and the location of the control point. 
First, the type of network must be specified, as follows: 

NTS(K) = 0 if the Kth network is source-free 

= 1 if sources are distributed on the Kth network as 
described in Section 1.2 

NTD(K) = 0 if the Kth network is doublet- free 

= 12 if doublets are distributed on the Kth network 
as described in Section 1.2. 

= 18 for a "type 18" wake network (see Section 1.3.2} 
= 20 for a "type 20" wake network (see Section 1.3.2) 

Both NTS and NTD must be specified for each network (K=l, ...., NNETT). If 
either parameter is zero, or if NTD = 18 or 20, one and only one boundary 
condition must be set at every panel center control point. Otherwise, two 








^T 
B i 



conditions are required at every panel -center control point. If NTD^O 
one and only one boundary condition is required at each network-edge 
control point, as shown in Fig. 1.10. 

The general form of the boundary condition. Equation (2) in Section 1.3, 
requires the specification of parameters CU, TD, DU, CL, IT, DL, and BET. 
This is accomplished by setting up a DO-loop in subroutine INPUT which 
assigns to every control point values of the following five parameters: 

ICA = A counter to distinguish one control point 

from another. 

NL0PT1, NL0PT2 = Parameters governing left sides for the two 

boundary conditions at the control point (if only 
one boundary condition is to be applied, one of 
these parameters is set equal to zero. 

NR0PT1, NR0PT2 = Parameters governing corresponding right sides of 

the two boundary conditions. 

That is, the left and right side of Equation (2) are controlled by setting 
values of two parameters, NLOPT and NROPT. This is done twice for each 
control point, the suffixes 1 and 2 being used to distinguish the two 
boundary conditions. 

The implications of particular choices of NLOPT and NROPT and NROPT are 
summarized in Tables 2.2 and 2.3, respectively. Options NLOPT = 2,3,4 and 
5 are used when normal mass-flux boundary conditions are to be applied, 
while NLOPT = 11, 12, 13 and 14 are useful for normal velocity boundary 
conditions. With any of these eight choices of NLOPT, the standard boundary 
condition of impermeability requires NROPT = 3 to select the right side 
of the boundary condition. NLOPT = 0 is used for a null boundary condition 
when only one is required. Thus, the Morino-type boundary conditions 
discussed in Section 1.3.1 require: 
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NL0PT1 = 5 
NROPT1 = 3 
NL0PT2 = 7 
NR0PT2 = 2 

at all panel-center control points, and 
NL0PT1 = 0 
NR0PT1 = arbitrary 
NL0PT2 = 7 
NR0PT2 = 2 

at all network- edge control points. For the wake networks discussed in 
Section 1.3.2, we set 
NL0PT1 = 0 
NR0PT1 = arbitrary 
NL0PT2 = 2 
NR0PT2 = 3 

at all control points (which are all on the network leading edge). 

For the superinclined networks of Section 1.3.3, 

NL0PT1 = 3 
NR0PT1 = 2 
NL0PT2 = 7 
NRQPT2 = 2 

at panel -center control points. At network-edge control points, 

NL0PT1 = 0 
NR0PT1 = arbitrary 
NL0PT2 = 7 
NR0PT2 = 2 

While the above a 1 ternati ves cover the most likely possibilities, the user 
may need to use other forms of the boundary condition (2). By setting 
NLOPT = 1, the coefficients of the left side may be set in subroutine 
INPUT. In so doing, one must also set the value of a parameter NCT for 
the control point in question: 
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If NLOPT = 1, set 

NCT =1 if DU = DL = TU = TL = 0 (i.e., iT the boundary condition 
will involve only the normal component of perturbation mass 
flux) 

= 4 if CU = CL = TU = TL = 0 (i.e., if the boundary condition 
will involve only the perturbation potential) 

= 2 otherwise 

Setting NL0PT1 = 1 requires the user to specify, in INPUT, the values of 
NCT1 , CU1 and CL1 (unless NCT1 = 4), DU1 and DL1 and (unless NCT1 = 1), and 
TU and TL1 (only if NCT! = 2). Setting NL0PT2 = 1 requires, instead, the 
specification of NCT2, CU2 and CL2, etc. 

Similarly, the right side of Equation (2) can be sec according to the user’s 
special needs in INPUT by setting NROPT = 1. If N0RPT1 = 1, the user 
specifies BET1, while if NR0PT2 = 1, it is BET2 that must be supplied in 
INPUT. To allow for the possibility that BET depends on the angles of 
attack and/or side-slip, or on the magnitude of the freestream velocity, 

BET! and BET2 are vectors of length 8. Thus, if NACASE 1, the user must 
set BET1(I) or BET2(I) for I = 1,..., NACASE. 

If the coefficients of the boundary condition (2) depend on the location 
of the control point, it will not be possible to set them in INPUT, since 
the locations of the control points are not set until after subroutine 
IBTRNS (see below) is called. In such cases, and possibly others, the 
user must supply subroutines (CCOF and/or CBET. These will be called by 
the program at time boundary conditions are being applied. 

Subroutine CCOF will have as arguments the parameters of the left side of 
equation (2), CU, CL, DU, DL, TU, TL as well as the indicator NBIN. The 
parameters TU and TT must be dimensioned variables, with dimension 3. 
Similarly CBET has arguments BET and NBIN where BET is a dimensioned 
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variable with dimension 8. Data on the control point location can be 
communicated through COMMON block /CNTRQ/. Should BET depend on the 
freestream flow conditions, they can be furnished CBET through COMM0N 
block /ACASE/. The structure of these blocks is given below. If both 
boundary conditions at a control point must be specified using CCOF 
and/or CBET they can be distinguished by the indicator NBIN=1 or 2. 

If NLOPT = 10, user supplies 

SUBROUTINE CCOF ( CU , CL, TU, TL, DU, DL, NCT) 

DIMENSION TU ( 3 ) ,TL( 3) 

If NROPT = 5, user supplies 
SUBROUTINE CBET (BET) 

DIMENSION BET (8) 

If either or both cases, the following common blocks may be userful: 
C0MM0N/CNTRQ/ZC(3) , ZNC{3), ZDC,IPC,ICC,JZC,JCN,KC 
C0MM0N/ACASE/ALPHA(8) , BETA(8, FSVM(8)) FSV{3,10),IACASE, NACASE 
C0MM0N/BC0ND/DUM (18), NCT, NLOPT, NROPT 


Here 


ZC(I) = Ith coordinate of control point whose boundary 

condition is being set 

ZNC(I) = Ith coordinate of outward surface normal at 

control 

KC = Network of control point 

FSVM(I) = Magnitude of freestream velocity in Ith case 

FSV(I,J) = Ith component of freestream velocity in Jth case 

and the other quantities have been described 
previously. 

and the other quantities have been defined previously. 
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If Morino-type boundary conditions are applied, the user gains one 
additional option, namely the method to be used for computing the flow 
velocity on the network surface. Generally, this can be done by summing 
the velocity fields due to the source and doublet distributions over all 
the networks. However, as noted in Section 1.3.1, use of Morino-type 
boundary conditions in theory makes the velocity on a particular panel 
identical with the gradient of the doublet strength on that panel. One 
can thus avoid calculation of the velocity influence coefficients, and 
so save substantial CPU and 10 time, or compute the surface flow velocities 
both from a sum over the entire system or networks and from the doublet- 
strength gradient, so as to gain a measure of the adequacy of the paneling 
The controlling parameter is IPOT(K), K being the network label. 

IPOT(K) * 2 if the velocities and pressures at control 

points on the upper surface of the Kth network 
are to be computed from the gradient of the 
local doublet strength. 

= 1 if they are to be computed both from the 
doublet strength gradient and by integrating 
over the entire system of singularities 

= 0 all velocities and pressures on the Kth network 
are to be computed only by integrating over the 
entire system of singularities 

= -1 if the velocities and pressures at control 
points on the lower surface of the Kth network 
are to be computed both by differenting the local 
doublet strength and integrating over the networks. 

= -2 if the lower surface velocities and pressures 
are to be computed by differentiating the local 
doublet strength. 



Options IPOT(K) = -1 and -2 would be used if it were the lower surface of 
the network which is wetted, and the potential on its upper surface set 
equal to zero (i.e., if NLOPT1 = 5, NR0PT1 = 6, NL0PT2 = 6, NR0PT2 = 2). 

If IPOT(K) = + 1, the forces and moments due to the upper or lower surface 
pressure distribution is calculated from the pressures corresponding to 
the doublet-strength gradient velocity field. 

Once the boundary conditions are set at a particular control point, they 
are filed by calling subroutine IBTRNS (ICA, CUI). As illustrated in 
Fig. 2.1, this is done on a networks-by-network basis with ICA being a 
cumulative control point counter. 

For non-wake networks, the control point index (IC in Fig. 2.1) ranges 
from 1 to (NN(K) + I) * (NM(K) + 1). For type 18 wake networks, it goes 
from 1 to (NN(K) + 1), while for type 20 networks IC goes from 1 to 1. 

2.5 EXAMPLES 

An example of subroutine INPUT is given in Fig. 2.2. 

This example was used to generate the data reported in Section 4.6 below. 
It could be adapted to any case for which the network geometries were 
available on punched cards in the following format: 

1. The first card gives the number of networks. 

2. For each network, the first card gives the numbers of rows 

and columns. For each column, succeeding cards give the global 
coordinates (x, y, z) of the grid points in order of increasing 
grid-point labels, two points to a card (the format for all 
cards is 6F10.6). This is done for each column in succession, 
so that the points are read in the order illustrated in Fig. 1.6. 
The example does assume Morino-type boundary conditions on all 
networks. Since the boundary-condition specification tends to 
be the most complex part of INPUT to the program, we give a 
number of examples below to illustrate what must be done in 
different situations. 
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SUBROUTINE INPUT 
C 

C THE COMMON BLOCKS WHICH FOLLOW MUST ALWAYS BE PART OF 
C THIS SUBROUTINE, VERBATIM 
C 

CQMM0N/ACASE/ALPHA(10) ,BETA( 8 ),FSVM( B )*,FSV{3, 8 ) , IACASE.NACASE 
COMMON/BCON/CU1 ,CL1 ,TU1 (3) ,T, 1 (3) ,DU1 ,0L1 ,BET1 ( 8 ) .NCT1 .NLOPT1 , 

+ NR0PTl,CU2,CL2.TU2(3),TL2(J),DU2 t DL2 I BET2( 8 ) ,NCT2 ,NL0PT2 ,NR0PT2 
COMMON/CASE/ ICASE.NCASE 

COMMON/COHPRS/ AMACH . BETAMS . fit r AM , SBETAH f ABETMS , ALPC t BETC , 

+ COMPD9 3 ) * AROTC ( 9 ) , AROTC 1(9) 

COWION/OATCHK/INDTCHK 

COMMON/ INDEX/NTS (50) . NTD ( 50 ) ,NH(50) ,NN(50) ,NZ(5C) »NP(50) ,NSS(50) , 

♦ NSD(SO) ,NC( 50) , N3C< 50) .NZAfSl ) ,NPA(51 ) ,N$$A(51 ) ,NSDA( 51 ) , 

♦ NCA(51 ) ,NBCA(51 ) , I POT (50) ,NNETT t NZMPT,NPANT,NSNGT,NSNGU,NSNGK, 

♦ NCTRT.NBCOT 
COMMON/MSPNT 5/ ZM( 3 , 1 500 ) 

COMMON/PRNT/ IGEOMP , ISINGP, ICONTP, IBCONP, 1EDGCP,ISINGS , IPRAIC 
COMMON/ SYMM/NSYMM 
C 

C THE FOLLOWING COWOH BLOCKS ARE OFTEN USEFUL 
C 

C IN /NCONS/ , PI-3.11. ,.,PI2-PI*2,PI4I-Pl/4 

C /SKRCH1/ PROVIDES SPACE FOR WHATEVER PREPROCESSING OF GRID POINT 
C DATA MAY BE REQUIRED 
C C0MM0N/NC0NS/PI.PI2.PI4I 

COMMON/SKRCHl /DUMMY (22500) 

C 

c SET TERMINATION LIMIT 

C NDTCHK - 0 

C 

C SET OUTPUT CONTROLS 
C ISINGP-0 $ IGEOMP-O 

1C0NTP-1 1 IBCONP-1 $ IEDGEP-1 $ ISINGS-1 S IPRA2C-1 

C 

C SET BASIC FLOW PARAMETERS 
C 

AMACH-2.05 
ALPC-2. S BETC-O, 

NACASE-1 

DO 100 IACASE-l.NACASE 
ALPHA(IACASE)-2. $ BETA(IACA5E)*0. 

100 CONTINUE 
NSYMM-1 
C 

C SET UP GRIO POINTS 
C 

C HERE THEY ARE READ IN IN STANDARD FORM USEO IN THE 
C 230 GRID GENERATION PROGRAM 


og 

§1 

Q <v 

El 

as 


C 

c 

C FIRST READ IN NUMBER OF NETWORKS 
C 

READ (5*1000) ANET 
1000 FORMAT(6F10.6) 

NNETT'ANH 

L-0 

DO 210 M.HNETT 


C 

C 

C 


c 

C 

c 


c 

C 

C 

C 

C 

c 


FOR EACH NETWORK, READ IN NUMBERS OF ROWS AND COLUMNS 
READ (5,1000) A,B 

H-A S N-B $ NM(K)-M $ NN(K)-N $ NPOINT-N*M 
DC 200 J-l ,N 

FOR EACH COLUMN, READ IN COORDINATES OF EACH GRID POINT 

IL-(J-1)*M*L+1 S IU-J*M*L 

READ (5,1000) ( (ZM( IC,I ) , IC-1 ,3) , I-IL,IU) 

200 CONTINUE 
210 L-NPOINT+L 

SET UP BOUNOARY CONDITIONS 

HERE THEY ARE PHI-0 ON LOWER SURFACE AT ALL CONTROL POINTS ANO 
SOURCE STRENGTH CANCELLING FREESTREAM AT ALL SOURCE CONTROL POINTS 


ICA-0 

DO 340 K-l ,NNETT 

NTS{K)«1 $ NTD(K)-12 $ IP0T(K)-1 

NMP1"NM(KM S NNP1-NN(K)*1 
DO 330 N-l.NNPl 
DO 320 M-1.NMP1 

NL0PT1-5 S NROPT.^3 $ NL0PT2-7 % NROPT2-2 

IF f (M.Fn. 1 1 .OR. (M.EQ.NHP1 ) .OR. (N.EQ.l ) .OR. (N.EQ.NNP1 )) NL0PT1- 


310 ICA-ICA+1 

CALL IBTRNS(ICA.CUl) 

320 CONTINUE 
330 CONTINUE 
340 CONTINUE 
C 

C THIS COMPLETES THE INPUT PHASE. NOW DIRECT THE PROGRAM TO CARRY OUT 
C THE CALCULATIONS ANO OUTPUT THE DESIREO RESULTS 
C 

RETURN 

ENO 



H 

r 

l 



Figure 2.2 Sample Input Subroutine 



2.5.1 Example 1: Morino-Type Boundary Conditions . 

The boundary-conditions section of the INPUT subroutine shown in Fig. 2.2 
assigns Morino-type boundary conditions to control points which are 
ordered column by column, just as are the grid points. As can be seen 
from Fig. 2.3, for a network with NM rows and NN columns, of grip points, 
there are NM + 1 rows and NN + 1 columns of control points. In the 
first and last rows and columns, only one boundary condition is required 
(5^= 0), while both boundary conditions are required at all other control 
points. This is reflected in the IF-statement of the programming 
illustrated in Fig. 2.2. 

For networks in which one edge collapses to a point, see Fig. 2.4, the 
number of control points is less than (NM + 1) ' (NN + 1). However, the 
programming listed in Fig. 2.2 will still be usable. The program recognizes 
the apparent redundancy in control points, and automatically stores the 
appropriate data for the tip control point just once. 

2.5.2 Example 2: Source Panel Network . 

The network consists only of sources, whose strength is to be specified 
equal to Q, a constant. It is still necessary to cycle the DO-loop thru 
the potential edge control points, so that the only changes from Example 1 
are that: 

NL0PT1 = 5 
NR0PT1 = 1 
NL0PT2 = 0 
NR0PT2 = 0 

and, before calling IBTRNS, one sets BET1 = Q. 

2.5.3 Example 3: Doublet Panel Network . 

The network consists only of doublets, whose strength is to be specified, 
and there is no adjacent network (such a network was specified to generate 
a disturbance in the checkout of the superinclined network). Some care 
must be taken in this specification. The doublet strength must vanish on 
any free network edge; otherwise a line vortex will be generated there. 



ORIGINAL PAGE IS 
OF POOR QUALITY 


NM(IO 


Figure 2.3 Control Points for Example 1: Quadrilateral 

Network 


Figure 2.4 Control Points for Example 1 : Triangular 


and the analysis specifically deletes line-vortex contributions to the 
velocity field (for numerical convenience). This is an example in which 
the right side of the boundary condition will depend on the control 
point location, and so a case where we must set NROPT = 5 and supply a 
subroutine CBET. 

Changes required in the DO-loop of Example 1 are as follows: 

DO 330 N = 1, NNP1 

DO 320 M = 1, NMP1 

NL0PT1 = 9 

NR0PT1 = 5 

NL0PT2 = 0 

NR0PT2 = 0 

ICA = ICA + 1 

CALL IBTRNS (ICA, CU1) 

320 CONTINUE 
330 CONTINUE 


Here is an example of the subroutine CBET which must be supplied: 
SUBROUTINE CBET (BET) 

DIMENSION BET (8) 

C0MM0N/CNTRQ/ZC(3), ZNC(3), DUM, IDUM(4), KC 
COMMON/DOUBLET/5IZE(2) 

DEL = SIZE (1) ** 2 
Y = ZC (2) 

Z = ZC(3) 

BET = 9./16. * (1. - (Y/DEL) ** 2) * (1. - (Z/DEL) ** 2)/DEL 

** 2 

RETURN 

END 


Note the communication of some miscellaneous data (the vector SIZE (2)) 
to CBET through yet another common block. 


ilia 


; n : 


2.5.4 Example 5: Multi -Network Wing Body Configuration . 

Our final example, the piece de resistance, is a wing-body problem. The 
configuration is shown in Fig. 2.5. As noted in Section 1.2.3, it is 
essential to make the points at which the leading and trailing edges of 
the wing meet the body at corner points of the singularity network. Thus, 
it is necessary to split the body surface into 6 networks, and the wing 
surface into two, as indicated in Fig. 2.6. Since this configuration is 
being to modeled for supersonic flow, it is not necessary to close off 
the base of the body to form a closed surface. In subsonic flow this 
would be necessary. The wing fillet, which is part of the body networks, 
gave us considerable trouble in insuring a perfect match between the grid 
points on the edges of the neighboring networks. Our solution is depicted 
in Fig. 2.7. 

It is necessary to model both the body's and the wing's wake. Because of 
the fillet, we required a "type 18" network to trail from the wing itself 
and a "type 20" to fill the gap between the wing's wake and the body. Two 
“type 18" networks were used to model the body's trailing wake. The 
complete system of networks is shown in Fig. 2.8 with the directions of 
the M and N vectors for each network given in the 3-views of Fig. 2.9. 

We are now ready to set up the boundary conditions. On the four wake 
networks (9, 10, 11, 12), the boundary condition is that the total mass 
flux be parallel to the surface. This is accomplished by setting NL0PT1 = 0, 
NRQPT1 = 0, NL0PT2 = 2, NR0PT2 = 3 for each control point on the wake 
network. Recall from Section 1.3.2 that networks only have control points 
at their leading edge. The configuration networks form a closed impermeable 
boundary so we wish the source strength to cancel the normal component of 
the freestream mass flux, and the doublet distribution to make the pertur- 
bation potential inside the configuration vanish. This is accomplished 
by setting NL0PT1 = 5, NR0PT1 = 3, NL0PT2 = 7 and NR0PT2 = 2 at all panel - 
control points ((W u - W & ) .IT = . "n and 0 A = 0) and setting NL0PT1 = 0, 

NL0PT2 = 7 and NR0PT2 = 2 at all panel edge control points (0^ = 0). 
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Aspect rd! io “ 1 .65 
iufjcr ratio * 0 10 
Angle of incidence “ 0 0" 

Reference arr*a “ 0 6256 (6 734 f i 

Span ” 101 .600 cm (40.00 in.) 

Body length, L » 227.866 cm (89.711 in ) 
Streamwise tip 


0 75 chord 


MS 194*310 
(76.500) 


MS 87.760(34.551) 


Norma) fence anil 
pitching moment 
gages -top and 
bottom — -r^ 


D 0.3937x * 1 .0096x^/1 


Constant diameter 8.748(3.444) 


Winy fLderpnee plane 


Centerline 


13.429 dia (5.287) 


3.149 

(1.240) 


' — Winy contour at BL 4.374 (1.722) 
Wing contour at centerline 


All dimensions in centimeters (inches) 


Arrow- Wing Body: General Arrangement 


FIGURE 2.5 
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The subroutine INPUT used for this configuration is illustrated in Fig. 2-10, 
In order to reduce computing costs, IPOT(K) = 2. This instructs the 
program to compute only the potential influence coefficients deleting 
the velocity influence coefficient calculations. With this option the 
velocities and pressures of control points are calculated by differenti- 
ating the local doublet strergth. 

The paneling coordinates have been generated external to the pilot code. 
Geometry coordinates are read in by network in 6F10.0 format. 

Networks numbered 10 or greater were inadvertently generated with the 
wrong NxM sense. A call to subroutine MNSWCH reverses the order of the 
input array. 

Networks 1 to 8 are designated composite source-doublet networks; networks 
9 and greater are designated wake type 18 networks, except for necwork 10 
which is a type 20. Morino-type boundary conditions are specified for 
the configuration surface networks (1-8). Once the boundary conditions 
have been specified for all networks, a call to subroutine IBTRNS stores 
the data on a random access file using ICA as an index. 


SUBROUTINE INPUT 

CGMHON/ACASE/ALPHA (8) ,BETA(B) ,FSUN(8) ,FSU(3,8) .IACASE, NACASF 
C0WON/8C0N/CU1 .Cll JU1 (3) ,TLl{3) *DU1 ,011 ,BET1(8) .NCT 1 ,NL0PT1 , 
CNRQPT1 ,CU2 ,CL2 JU2(3) , TL2(3) ,DU2 ,DL2 .BET2(8) ,NCT2 ,NLOPT2 .NROPT2 
COHMON/CASE/ ICASE ,NCASE 

CQMMQN/COHPRS/AMACH.BETAMS ,BETAM,SBETAM .ABETMS ,ALPC,BETC, 
CC0MPD(3),AR0TC(9),AR0TCI(9),CZINU(3,3) 

COtmi/DATCHK/NDTCMK 

COMMON/ IND£X/NTS(5Q) ,NTD(50) ,NM{50) ,NN(50) ,NZ(50) ,NP(50) ,NSS(50) „ 
CNSD(50) ,NC(50JNBC( 50) ,NZA(51 1 , NPA(Sl) ,NSSA(5I ) ,NSDA(51 ) ,NCA(51 ) . 
DtBCA{51 ) .IPOT(SO) ,NNETT .NZHPT ,NPANT ,NSNGT,NSNGU,NSNGK,NCTRT .NBCOT 
CQMNQN/MSPNTS/ZM( 3,1530) 

COrtiON/NCONS/P I ,P 1 2 r P 1 4 1 

COftON/PRNT / IGEOMP , IS I N6P , ICONTP , IEDGEP , ISINGS, IPRAIC 
CGH10N/5YHN/NSYW 

C 

C SET OUTPUT CONTROLS 
C 

ISINGP-O 

IGCOMP-O 

e 

C SET-TERMINATION LIMIT 
C 


C 

C SET INDICATOR FOR COMPUTING POTENTIAL INFLUENCE COEFFICIENTS 
C ONLY TO REDUCE CPU AND 10 COSTS 
C 

IP0T(K)*2 

C 

C FOR EACH NETWORK, READ IN NUMBERS OF ROWS AND COLUMNS 
C 

R£A0(5,100)A,B 
100 FORMAT (6F10.0) 

105 H-A 
N-B 

NN(K)*N 

NM(k)-M 

C 

C READ NETWORK GEOMETRY 

C 

DO 20 J-1.N 

READ(5.100){(ZH(JJ,i+L),JJ-l,3M-\.N) 

L-L^M 

20 CONTINUE 

NZAPl-NZA(k)*l 

C 

C N X M ORDER OF NETWORKS 10 OR GREATER IS INCORRECT 
C CALL SUBROUTINE MNSWCH TO REVERSE ORDER 

C 

IFOC.GE. 10) CALL' MNSWCH{ZM{ 1 .NZAP1 ) ,NM{K) ,NN(K) ) 


NDTCHK-0 
C 

C sn BASIC FLOW PARAMETERS 
C 

AMACH-1.7 
NSYMM-1 
NACASE-1 
ALPHA(l)-2.0 
ALPC-2.0 
IETC-0 
C 

C SET UP NETWORK GEOMETRY 
C 

C FIRST REAO IN NUMBER OF NETWORKS 
C 

READ(5,100) ANET 
NNETT-ANET 
L-0 
ICA-0 

C 

C NZA NEED NOT BE SET IN INPUT. BUT IS USED HERE ONLY TO 
C FACILITATE SUBSEQUENT CALL Of HNSWCH 
C 

NZAOJ-0 

DO 3000 K*1 ,NNETT 
C 

C SET SINGULARITY TYPE 
C 
C 

C COMPOSITE SOURCE - DOUBLET NETWORKS 
C 

NTS(K)-1 
NTD(K)-12 
C 

C NETWORKS 9 or GREATER ARE WAKE NETWORKS 
C 

IF(K.GE. 9) NTS ( K ) *0 
IF(K,GE.9) NTD(K)-18 
IF(K.EQ.IO) NT0(K)*20 

Figure 2* 10 INPUT Subroutine for Arrow-Wing Body 
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C NZ NEED NOT BE SET IN INPUT. BUT IS CALCULATED HERE ONLY 
C TO FACILITATE SUBSEQUENT CALLS OF MNSWCH 
C 

NZ(K)-NM{K)*NN(K) 

NMPl ■NM( K ) ♦ 1 
NNPI-NN(K)*1 
C 

C SET BOUNDARY CONDITIONS 

c 

DO 150 N-l ,NNP1 
00 140. H-lp NMPl 
C 

C WAKE NETWORKS REQUIRE BOUNDARY CONDITIONS ONLY FOR LEADING EOGE 
C 

IF((NTD(K).EQ.l*}.AND.iM.GT.1)) GO TO 150 
IF((nTD(K).EQ.20).AND,(N.GT.1)) TO TO 3000 

C 

C COMPOSITE PANEL MASS FLUX BOUNDARY CONDITIONS 
C 

NL0PT1-5 

NR0PT1-3 

NL0PT2*7 

NR0PT2-2 

C 

C ONLY ONE BOUNDARY CONDITION SHOULD BE ENTEREO FOR EDGE CONTROL POINTS 
C 

IF((H.EQJ).OR.(M.EQ.NHPl).QR.(N.Eq.1).OMN.EQ.NNPl))NLOPTl-0 
130 ICA-lCA+1 
C 

C IBTRNS SETS UP BOUNDARY CONDITION DATA ON A RANDOM 
C FILE USING ICA AS AN INDEX 
C 

CALL IBTRNS ( I CA.CUl) 

140 CONTINUE 
ISO CONTINUE 

3000 NZA(K*1 )"NZA(K)+N2(K) 

RETURN 

END 
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QUANTITY 

DEFINITION 

1 DEFAULT 
i VALUE 

COMMON 

BLOCK 

ALPC 

Compressibility direction angle (of attack) in 
degrees 

0 . 

/C0MPRS/ 

ALPHA (I) 

Ith angle of attack in degrees 

0 . 

/ACASE/ 

AMACH 

Mach Number 

0 . 

/C0MPRS 

BETA (I) 

Ith angle of yaw in degrees 

0 . 

/ACASE/ 

BETC 

Compressibility direction angle (of yaw) in degrees 

0 . 

/C0MPRS/ 

BET1(I) 

Ith right hand side of first boundary condition 
(see BET of equation{2)) 

0 . 

/BC0N/ 

BET2 ( I ) 

Ith right hand side of second boundary condition 
(See BET of equation (2)} 

0 . 

/BC0N/ 

BREF 

x-axis value of T R (SECTION 1.4) 

1 . 

/FMC0F/ 

CL1 

Value of CL for first boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

CL2 

Value of CL for second boundary condition 
(See equation (2)) 

0 . 

BC0N/ 

CREF 

y-axis value of T R (SECTION 1.4) 

1 . 

/FMC0F/ 

CU1 

Value of CU for first boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

CU2 

Value of CU for second boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

DL1 

Value of DL for first boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

DL2 

Value of DL for second boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

DREF 

z-axis value of T R (SECTION 1,4) 

1 . 

/FMC0F/ 

DU1 

Value of DU for first boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

DU2 

Value of DU for second boundary condition 
(See equation (2)) 

0 . 

/BC0N/ 

FSVM(I) 

Ith Freestrea- velocity magnitude 

1 . 

/ACASE/ 

IBC0NP 

=1 for printout of boundary condition diagnostic data 

1 

/PRNT/ 

IC0NTP 

=1 for printout of control point diagnostic data 

1 

/PRNT/ 










QUANTITY 


DEFINITION 


DEFAULT 

VALUE 


COMMON 

BLOCK 


=1 for printout of edge matching diagnostic data 
“1 for printout of geometry diagnostic data 

* r 0 for computation of resultant velocities at 
i network control points using velocity 

influence coefficients 

±1 for added computation of resultant velocities 
at network control points using doublet 
< strength gradient 

±2 for computation of resultant velocities from 
doublet strength gradient 

(+ for upper surface velocities with lower 
surface perturbation stagnation and - for 
„ the opposite) 

= L for printout of influence coefficient diagnostic 
data for L^* control point 

=1 for printout of influence coefficient diagnostic 
data 

=1 for printout of singularity spline diagnostic dat 

=1 for printout of resultant values of singularity 
strength and gradient at panel corners, center and 
edge midpoints 

Number of simultaneous solutions (freestream vectors 

r 

= 1 first boundary condition concerns normal 

mass flux only 

' 2 first boundary condition is arbitrary 

4 first boundary condition concerns potential 
w only 

Same as NCT1 but for second boundary condition 

* JO program run to completion 

I 1 program stop immediately before influence 
^ coefficient calculation 

First boundary condition left hand side option 

Second boundary condition left hand side option 

First boundary condition right hand side option 

Second boundry condition right hand side option 




/PRNT/ 

/PRNT/ 


/INDEX/ 


/PRNT/ 


/PRNT/ 

/PkNT/ 


/ACASE/ 

/BC0N/ 



/BC0N/ 

/DATCHK/ 


/BC0N/ 

/BC0N/ 

/BC0N/ 

/BC0N/ 
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QUANTITY 


NM(K) 

NN(K) 

NNETT 


NSYMM* 


NTS(K) 

NTD(K) 


ZM(I,L) 


DEFINITION 


Number of grid point rows in Kth network 
Number of grid point columns in Kth network 
Total number of networks 


= 0 for no planes of flow symmetry 

1 for x-z plane of flow symmetry with y 0 
< half of configuration input 

2 for x-z and x-y planes of flow symmetry 

| with y ^ 0 and z .> 0 quarter of configuration 
k input 

Source type of Kth network 

Doublet type of Kth network 

Configuration reference area (SECTION 1.4) 

x-component of R R (SECTION 1.4) 

y-component of R R (SECTION 1.4) 

Ith component of Lth grid point (L is a cumulative 
index for all grid points in the total configuration) 

z-component of R R (SECTION 1.4) 


DEFAULT COMMON 
VALUE BLOCK 


/INDEX/ 






/INDEX/ 

/INDEX/ 


/SYMM/ 


/INDEX/ 

/INDEX/ 

/FMC0F/ 

/FMC0F/ 

/FMC0F/ 

/MSPNTS/ 


/FMCOF/ 


*If NSYMN * 1 or 2 is selected, the force and 
moment coefficients will be computed from the 
pressures acting on half, and a quarter of 
the configuration, respectively. 
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NLOPT 

Left Side of Equation 

Significance 

0 


Boundary condition is 
ignored (used when none or 
one boundary condition 
is to be applied to a 
control point 

1 

Arbitrary 

Parameters of boundary 
bondition are set in 
INPUT 

2 

w • n 

Normal component of 


u 

perturbation mass flux 
on upper surface 

3 

"a ' " 

Normal component of 
perturbation mass flux 
on lower surface 

4 

i(w n + w ) • n 

Average of upper and lower 


C U X* 

surface normal perturbation 
mass flux 

5 

(w M - wj • n’ 

Difference between upper and 


u 

lower-surface normal mass flux 
(source strength 

6 

4> 

Upper-surface perturbation po- 
tential 

7 


Lower-surface perturbation 
potentials 

8 

■k4>.. + <m 

Average of upper and lower- 


c U X, 

surface perturbation potentials 

9 

♦u “ 

Difference between upper and 
lower surface potential (doublet 
strength) 

10 

Arbitrary 

Parameters of boundary conditions 
reset in subroutine CCOF 

11 

v„ 

Upper -surf ace perturbation 


u 

velocity 



Eg 
















Table 2.2 (Cont.) 

Effect of Various Choices for KLOPT on 
Left side of boundary Condition (2} 


HLOPT 

Left Side of Equation 

Significance 

12 

V *. 

Lower-surf ace perturbation 
velocity 

13 

|( v u + v z> 

Average of upper and lower- 
surface perturbation 
velocities 

14 

< u u - v l> 

Difference between upper- 
and lower-surface perturba- 
tion velbcities 









1 - ±1 


1 . 


Table 2.3 

Effect of Various Choices for NROPT on 
Right Side of Boundary Condition (2) 


EROPT 


BET (Right Side of 
Equation (2) ) 

! 

Significance 

Arbitrary 

0 

- V • n 

co 

BET will be specified in INPUT 

If NL0PT=2 , equation (2) says 
total mass flux is parallel to 
surface. If NLOPT-5, equation 
(2) specifies source strength 
so as to cancel freestream 
mass flux 

-V c • (x/B 2 ,y,z) 

If KIiOPT=7 and NTS=0 (no sources) 

V c =freestream velocity 

it can be shown that equation 

in terms of components 

(2) implies no flow thru sur- 

along compressibility 
axes) 

face 

Arbitrary 

BET will be specified in CBET. 

V • n 

Same as 3 but for surface with 

oo 

lower side exposed to flow 










3.0 Output Data Description 


In this section we describe the organization of the output. The 
nomenclature employed in the output is summarized completely in Table 3-1. 

3.1 GEOMETRY (MESH POINT) DATA 

The program always lists the coordinates of all grid points, along with 
their unique labels, row number, column number, and network number. 
Specifying the network geometry is perhaps the trickiest part of setting 
up INPUT. The geometry data are always worth examination to see if any 
mistakes have been made, particularly in assigning grid points to rows 
and columns. These data are illustrated in Fig. 3.1. This is followed 
immediately by an identification of whatever subpanels are "superinclined" 
i.e., whose normals are inclined to the compressibility reference axis 
by more than the Mach angle. Except in special cases, such as represen- 
tative of inlets, superinclination of a panel indicates that the small- 
disturbance assumptions underlying the analysis are being violated. 
Moreover, the answer will be invalid if boundary conditions at the control 
point(s) on the panel assume subinclination. In addition, the program 
will identify subpanels that are close to being superincl ined (probably 
violating small disturbance assumption and print the value of the dot 
product of the subpanel normal with the subpanel conormal . 

The listing of the actual corner points is followed by a printout of 
which sides of which doublet networks the program believes abut. This 
data is shown in Fig. 3.2 and is particularly useful in helping the user 
to determine whether the input grid points satisfy the second point of 
Section 1.1.2. The data is grouped according to abutments (i.e., groups 
of coincident network sides). Network sides are numbered counter clock- 
wise from 1 to 4 with the first side corresponding to the first row of 
gria points and the fourth side corresponding to the first column of grid 
points (see Fig. 3.3). Of particular interest are abutments involving 
only one side whose control points are characterized as forcing doublet 




^ strength to vanish. This implies that the side abute nothing - which 

may possibly be contrary to the user‘s wish. (The characterization of 
control points does not apply to network control points.) 


3.2 CONTROL POINT DATA 

Unless the user specifies ICONTP ? 1, the program next lists, for each 
potential* control point, its label (what is called ICA in the INPUT), 
network number, panel number, subpanel number, coordinates, the components 
of the unit normal vector to the subpanel on which it is situated, and 
an edge matching indicator. The unit normal should be directed from 
lower to upper surface of the network. This is worth checking. If it 
is misdirected, a call to MNSWCH may be necessary; see Section 2.3. 
Excessive variation from one control point to the next of the components 
of n often indicates that the panel size is too large. These data are 
illustrated in Fig. 3-4. 

3.3 BOUNDARY CONDITION DEFINING PARAMETERS 

Unless the user specified IBCONP f 1, the program outputs, for every 
control point, the parameters CU, CL, TU, TL, DU, DL, NCT and BET for the 
boundary condition(s) imposed at that point. Also included are the 
control point label ICA and the values of NLOPT and NROPT. Especially in 
multi-network problems, for which the setting of boundary conditions on 
network edges can be tricky, these data are often worth a look. These 
data are illustrated in Fig. 3-5. 


3.4 PROBLEM AND NETWORK INDICES 

These data, illustrated in Fig. 3-6, are always generated. They give the 
basic statistics on the problem size, both totally and network by network. 
Current limits are: 


♦Dataware also output for the midpoints of panel edges which lie 
on the boundaries even of source-only networks; these points are 
not otherwise treated as control points. 
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Total number of grid points (NZMPT) ^ 1500 
Total number of control points (NCTRT) « 1500 
Total number of panels « 1000 

3.5 EDGE DOWNWASH CONDITIONS 

These may be suppressed by stipulating IEDGEP f 1. They give the calcul- 
ated matching information at control points on the edges of doublet 
networks. The output is intended to verify that this matching has 
occurred. These data are also shown in Fig. 3.6. 1 

3.6 EQUATION SOLUTION DATA 

These data, shown in Fig. 3.7 and Fig. 3.8 are always printed out. They 
give the order of the algebraic system governing the singularity strength, 
an error estimate, and pivoting information. 

3.7 AERODYNAMIC INPUT DATA 

Next the program always outputs the values of M, a , and 8 also shown 
in Fig. 3.8. 

3.8 SINGULARITY DATA 

Unless ISINGS is set f 1, the program next prints out the source strength, 
the doublet strength, and the components of the gradient of the doublet 
strength along the row and column directions and in direction normal 
thereto as well as along the global coordinate axes. This is done for 
nine points on every panel: the four corners, the mid points of each 

side, and the panel center. These data, shown in Fig. 3.9, are useful 
in checking the continuity of the doublet strength from one network to 
another. ' 

3.9 LOCAL AERODYNAMIC DATA 

These data are always printed out, and usually contain the major results 
of the analysis. For every panel center-control point on every network, 
the following quantities are listed as illustrated in Fig. 3.10: 
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Source Strength 
Doublet strength 
Doublet strength gradient 
Perturbation velocity potential* 

Total velocity potential* 

Total mass flux potential* 

Perturbation mass flux vector* 

Total mass flux vector* 

Normal and tangential components of total 
Mass flux vector* 

Normal component of perturbation mass flux vector* 

Pressure coefficient* 

Quantities starred are listed for both the upper and lower surfaces of 
the network. Four different pressure formulas are used for the pressure 
coefficient: 


C p = - 2u (Linearized) 

C - - 2u - v 2 - w 2 (Slender Body) 

r 

Cp = - 2u - (1 - Lf ) u 2 - v 2 - w 2 (Second Order) 

C = {[1 + M 2 (l - u 2 - v 2 ~ w 2 )]* - 1} (Isentropic) 

P y H 2 ^ 

CO 



Here ( u, v, w ) is the perturbation velocity vector, referred 

to compressibility coordinates. Also listed are the differences in the 
four pressure coefficients across the network (upper surface value minus 
lower), the control point label, the panel label, the control point 
coordinates, and the components of the freestream velocity. If the 
INPUT parameters IPOT(K) = 1 for the Kth network velocity potentials, mass 
flux vectors, and pressure coefficients on the upper surface are computed 
(and listed) twice. The first time (line 2 or the output block), their 
values are based on the local gradient of the doublet strength; the second 
time (line 3 of the output block), on the integral over the entire 
singularity distribution. If IP0T(K)= -1 the same holds for lower surface 
quantities with the lower line being that data computed from the doublet 
strength gradient. For the options IPOT(K) = + 1 there should be little 
difference between the redundant data. If a significant difference does 
occur it usually means that paneling is inadequately dense or the surface 
is insufficiently closed (i.e., there is too large a hole in the paneling 
through which the fluid is flowing) for the use of the Morino-type 
boundary conditions. An almost certain indicator of such inadequacies 
is the existence of significant non-zero normal flows (i.e., WNU or WNL 
as computed from the influence coefficient) on the side of the surface 
bounding the physical flow of interest. 

These results are followed by listings, illustrated in Fig. 3.11, of 
the forces and moments acting on the network surface both totally and 
for each column of the network. The three lines of data given for each 
column (and for the whole network) refer to the forces and moments on 
the upper surface, the lower surface, and the sum of the two. The user 
may have to combine contributions from different networks manually rather 
than reiy on the automatically accumulated totals. For example, if the 
wing of a wing-body combination has no thickness, the wing’s contribution 
to the force involves botii its upper and lower surfaces, while only the 
upper (outside) surface of the body makes a physically meaningful 
contribution. 
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FIGURE 3.5 BOUNDARY CONDITION DEFINING PARAMETERS 
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FIGURE 3.6 PROBLEM AND NETWORK INDICIES AND EDGE DOWNWASH CONDITIONS 
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FIGURE 3.9 SINGULARITY GRID DATA 
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TABLE 3-1 DEFINITION OF OUTPUT QUANTITIES 



DEFINITION 

WHERE 

WHEN 

QUANTITY 

PRINTED 

PRINTED 


Problem and Network Indices 



NNETT 

Total number of networks 

R0W 

ALWAYS 

NZMPT 

Total number of grid points 

FL0W 

ALWAYS 

NPANT 

Total number 'of panels ORIGINAL PAGE IS 

FL0W 

ALWAYS 

OF POOR QUALITY 

FL0W 

ALWAYS 

NSNGT 

Total number of singularity parameters 

NSNGU 

Total number of unknown singularity parameters 

FL0W 

ALWAYS 

NSNGK 

Total number of a-priori known singularity 
parameters 

FL0W 

ALWAYS 

NCTRT 

Total number of control points 

FL0W 

ALWAYS 

KBC0T 

Total number of boundary conditions 

FL0W 

ALWAYS 

NTS 

Network source type 

FL0W 

ALWAYS 

NTD 

Network doublet type 

FL0W 

ALWAYS 

NK 

Number of corner point rows in network 

FL0W 

ALWAYS 

NN 

Numbeh of corner point columns in network 

FL0W 

ALWAYS 

NZ 

number of corner points in network 

FL0W 

ALWAYS 

NP 

Number of panels in network 

FL0W 

ALWAYS 

NSS 

Number of singularity parameters of network 
source spline 

FL0W 

ALWAYS 

NSD 

Number of singularity parameters of network 
doublet spline 

FL0W 

ALWAYS 

• NC 

Number of control points in network 

FL0W 

. ALWAYS 

NBC 

Number of non- null boundary conditions in 
network 

FL0W 

ALWAYS 


. Near Field/Far Field Counters 



NFS 

Number of near field source AIC calculations 

FL0W 

ALWAYS 

NFD 

Number of near field doublet AIC calculations 

FL0W 

ALWAYS 

FFS 

Number of far field source AIC calculations 

FL0W 

ALWAYS 

FFD 

Number of far field doublet AIC calculations 

FL0W 

ALWAYS 
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TABLE 3-1 DEFINITION OF OUTPUT QUANTITIES (CONT'D) 
DEFINITION 

Equation -Solution Data 


This data is self explanatory. The printout 
begins with the words "BEGIN S0LUTI0N 0F 
mm = m and ends with the "S0LUTI0N 
C0MPLETE" 


Singularity Data 
Solution singularity parameters 


Cumulative index of panel on which singularity 
distribution is evaluated 


Local row index uf evaluation point 
Local column index of evaluation point 
Global coordinates of evaluation point 
Source strength value at evaluation point 
Doublet strength value at evaluation point 


Global coordinates of gradient of doublet 
strength at evaluation point 


Derivative of doublet strength in (row, column) 
directions respectively 


Derivative of doublet strength in directions 
normal to (row, column) directions respectively 


Local Aerodynamic Data 

Cumulative control point index 

' « * 

Index of panel containing control point 

Global coordinates of control point 

Doublet strength 

Global coordinates of doublet gradient 
Source strength 

Freestream velocity vector in global coordinates 


WHERE WHEN 
PRINTED PRINTED 


ALWAYS 


0UTPUT ALWAYS 
0UTPUT ISINGS=1 


0UTPUT ISINGS=1 
0UTPUT ISINGS=1 


0UTPUT ISINGS=1 
0UTPUT ISINGSrl 


0UTPUT ISINGS=1 
0UTPUT ISINGS=1 


0UTPUT ISINGS=1 
0UTPUT ISINGS=1 


0UTPUT 

0UTPUT 

0UTPUT 

0UTPUT 

0UTPUT 

0UTPUT 

0UTPUT 


ALWAYS 

ALWAYS 

ALWAYS 

ALWAYS 

ALWAYS 

ALWAYS 

ALWAYS 
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QUANTITY 


(WXU.WYU, 

WZU) 


(PWXU,PWYU, 

PWZU) 

CPLINU 

CPSLNU 

CP3NDU 

CPISNU 


(WXL,WYL, 

WZL) 


(PWXL»PWYL, 

PWZL) 

CPLINL 
CPSLNL 
• CP2NDL 
CPISNL 


PHIUI 
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TABLE 3-1 DEFINITION OF OUTPUT QUANTITIES (CONT'D) 


DEFINITION • 

Upper surface total potential 

Upper surface total mass flux vector in global 
coordinates 

Upper surface perturbation potential 

Upper surface perturbation mass flux vector in 
global coordinates 

Upper surface linearized pressure coefficient 
Upper surface slender body pressure coefficient 
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WHERE I WHEN 
PRINTED PRINTED 

0UTPUT i ALWAYS 
0UTPUT ALWAYS 


{3UTPUT ALWAYS 
0UTPUT ALWAYS 


0UTPUT ALWAYS 
0UTPUT ALWAYS 


Upper surface second order pressure coefficient 0UTPUT 1 ALWAYS 


Upper surface isentropic pressure coefficient 

Lower surface total potential 

Lower surface total mass flux vector ip. global 
coordi nates 

Lower surface perturbation potential 

Lower surface perturbation mass flux vector in 
global coordinates 

Lower surface linearized pressure coefficient 

Lower surface slender body pressure coefficient 

Lower surface second order pressure coefficient 

Lower surface isentropic pressure coefficient 

Normal component of upper surface total mass 
flux vector 

Normal component of lower surface total mass 
flux vector 

Magnitude of tangential component of upper 
surface total mass flux vector 

Magnitude of tangential component of lower 
surface total mass flux vector 

Upper surface total mass flux potential 

Lower surface total mass flux potential 

Normal component of upper surface perturbation 
mass flux vector 


0UTPUT ALWAYS 
0UTPUT ALWAYS 
0UTPUT ALWAYS 


0UTPUT ALWAYS 
0UTPUT ALWAYS 


0UTPUT ALWAYS 
0UTPUT ALWAYS 
0UTPUT ALWAYS 
0UTPUT ALWAYS 
0UTPUT ALWAYS 


0UTPUT ALWAYS 


0UTPUT ALWAYS 


0UTPUT ALWAYS 


0UTPUT ALWAYS 
0UTPUT ALWAYS 
0UTPUT ALWAYS 
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QUANTITY 

DEFINITION 

WHERE 

PRINTED 

WHEN 

PRINTED 

. PHNL 

Norms! component of lower surface perturbation 
mass flux vector 

0UTPUT 

ALWAYS 

CPLIND 

Difference between upper and lower surface 
linearized pressure coefficient 

0UTPUT 

ALWAYS 

CPSLND 

Difference between upper and lower surface ' ~ 
slender body pressure coefficient 

0UTPUT 

ALWAYS 

CP2NDD 

Difference between upper and lower surface • 
second order pressure coefficient 

0UTPUT 

ALWAYS 

CPISND 

Difference between upper and lower surface 
isentropic pressure coefficient 

0UTPUT 

ALWAYS 


Force and Moment Data 



AREA 

Total area of panels 

FKCAL 

ALWAYS 

. (FX„FY,FZ) 

Global coordinates of force coefficient 
{Upper surface, lower surface, difference) 

FMCAL 

ALWAYS 

(MX.HY.MZ) 

Moment coefficients about global principal 
axes (upper surface, lower surface, difference) 

FMCAL 

ALWAYS 
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4.0 Computational Results 

In this section we present results obtained by applying the advanced 
panel pilot code to a number of different problems chosen to test its 
range of applicability and to obtain advice on how to use it most 
effectively. In most of these cases, exact solutions, or at least 
analytical results based on linear theory, are available for purposes of 
comparison. For the more complex cases comparisons have been made with 
experimental results. 

4.1 AXISYMMETRIC FLOW PAST CIRCULAR CONES 

One of the simplest cases of supersonic flow for which an exact solution 
is available is the flow past a cone at zero angle of attack. For this 
case the surface pressure is constant and depends only on the upstream 
Mach number M m and cone angle 0 c - Results obtained with the panel 
code are compared with the exact solution in Fig. 4.1. 

Both velocity and mass flux boundary conditions were tried. Requiring 
( w + ) . n = 0 (mass flux boundary condition) on the boundary is 

seen to give much better results than if (v + } . n = 0 (velocity 

boundary condition) on the surface when the Mach number and cone 
angle 0 C are relatively small; i.e., within the limits we expect for a 
linearized theory. However, the results deteriorate less rapidly with 
increasing Mach number when velocity boundary conditions are used. 

The simplicity of flow past a cone make it a suitable vehicle for a 
study of convergence. Fig. 4.2 shows the effect of the number of 
meridianal panels employed on the cone surface pressure at three Mach 
numbers. In every case three rows of panels were used in the axial 
direction (very little variation of pressure was observed in that 
direction). 
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FIGURE 4.1 PRESSURE ON CONE AT 0° ANGLE OF ATTACK 










It can be seen that the results are essentially converged when the angle 
subtended by the panels is less than 10°. This holds even when the results 
are converging to the incorrect answer because the Mach number and/or cone 
angle is too large. 

Somewhat faster convergence was obtained when the grid points were selected 
outside the cone surface so that the panels "interpolated' 1 rather than were 
inscribed by the cone, see Fig. 4.3. While such a device is not generally 
practical in using the panel method, we did take advantage of it in our 
studies related to cones. This modeling yields better results because the 
cross-sectional area of the "interpolated" cone better matches the actual 
cross-section area than does the area of an inscribed cone. 

4.2 FLOW PAST CIRCULAR CONES AT ANGLE OF ATTACK 

Extensive tabulations are now available of the exact solution of the flow 
past cones at angle of attack, Reference 3. As would be expected, the per- 
formance of the panel code deteriorates with increasing angle of attack for 
a given Mach number and cone angle. Fig. 4.4 thru 4.7 show the effects of 
the choices of compressibility axis and pressure formula for two conditions, 
one within the range of accuracy of the linearized theory (M m = 1.5, J3 C = 10°, 
a = 5°) and one beyond it (M m = 2.0, 0 C = 10°, a = 10°). These data 
clearly show the desireabil ity of aligning the compressibility axes with the 
freestream direction (i.e., setting « = a c 0 = 8 C ) and indicate that 
either the isentropic or slender-body pressure formulas are superior to the 
full second-order pressure formula. See Section 3.9 for the explicit formulas. 

Figs. 4.8 and 4.9 give another comparison of the use of velocity and mass 
flux boundary conditions. As was observed in connection with Fig. 4.1, mass 
flux boundary conditions give much better results when linear theory is ex- 
pected to be valid, but velocity boundary conditions are less wrong outside 
those limits of validity. 
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4.3 FLAT DELTA WINGS 

As indicated in Section 1.1, it is impossible to panel a triangular wing 
with exclusively quadrilateral panels, and one must choose the vertex in 
whose vicinity the panels will be triangular. Fig. 4.10 shows the effect 
of paneling on the velocity ratio Vy/V x on the surface of a flat delta 
wing with straight trailing edge at angle of attack. According to linear 
theory reference 4, this ratio ought to be linear in y/x. Fig. 4.10 
shows that the conical paneling "A" gives results in excellent agreement 
with linear theory, as would be hoped, since this paneling best accommodates 
the nature of the solution. However, paneling "C", in which the tri- 
angular panels are located near the midpoint of the trailing edge, works 
just as well, which makes the failure of paneling "B", in which the 
tips are treated with triangular panels, all the more surprising. This 
effect of paneling needs further study. 

4.4 DELTA WINGS WITH THICKNESS 

The flow about thick delta wings was also calculated. Here the upper 
and lower wing surfaces were modeled by separate composite panel net- 
works. The boundary conditions employed to obtain these results were 
the usual impermeable surface conditions described in Section 1.3 and 
2.5.1. 

Results for a case with a supersonic leading edge are shown in Fig. 4.11. 
These results are in good agreement with the conical flow linear theory. 

The computed results do show a slight pressure oscillation near the Mach 
cone from the leading edge. This is quite typical in that the panel 
method assumes a degree of continuity which is often absent in super- 
sonic flow. Fig. 4.12 shows the pressure distribution on a delta wing 
with a subsonic leading edge. Results are in excellent agreement with 
linearized theory. 
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FIGURE 4.11 PRESSURE DISTRIBUTION NEAR TRAILING EDGE OF 
DELTA WING WITH SHARP SUPERSONIC LEADING EDGE 
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The results plotted in Fig. 4.11 and 4.12 were obtained with a conical 
paneling like that called “A" in Fig. 4.10. Similar results were 
obtained when the triangular panels were situated near the midpoint of 
the trailing edge, as in paneling "C" of Fig. 4.10. 

4.5 FLOW INDUCED ON A PLANE BY A CONE ABOVE THE PLANE 
In order to assess the utility of the panel code to study interference 
effects, we considered the problem illustrated in Fig. 4.13, in which 
a cone at zero angle of attack is positioned above a plane, also at 
zero angle of attack. An exact solution can be constructed by the method 
of images, i.e., by taking the x-y plane to be a plane of symmetry. 

Quite good results were obtained for the pressure distribution on the 
plane when it was taken to be a source network with just one row of 
panels in front of the trace of the Mach cone on the plane. When doublets 
were used to meet the flow tangency condition on the plane, the results 
were much poorer, with noticeable "leakage" of the pressure distribution 
ahead of the Mach cone trace. This is due to the larger data base used 
to generate the doublet spline. As described in Appendix A, the doublet 
distribution over a panel is related to the doublet strength at 20 
surrounding panels, while the source strength spline relates the source 
strength at just 8 contigious panels. 

The errors shown in Fig. 4.13 for the pure-doublet network were notice- 
ably reduced (but still worse than for the purs-source network) both 
when more panels were put in ahead of the Mach cone trace and when, by 
refining the mesh spacing, more control points in the second row of 
panels (just aft of x=1.0, the most forward position of the Mach cone 
trace) were behind the Mach cone. These results are also consistent 
with the dependence of the doublet spline on the doublet strength on a 
fairly large number of contigious panels. 
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4.6 FLOW PAST A CAMBERED, TWISTED WING 

Extensive experimental data are available for the Carlson wing 2T 
reference 6, an arrow wing with sharp leading and trailing edges, twist, 
and camber. Results of the panel code are compared in Fig. 4.14 and 4.15 
with data taken at Mach 2.05, at which the leading edge is subsonic but 
the trailing edge supersonic. The theoretical results compare well with 
the experimental data except at the wing tip, where the theory tends to 
over predict the experimental results. This may be due to the paneling 
in this region, which left an open gap between upper and lower surfaces 
of the wing at the tip. The boundary condition imposed that 0^ = 0 
requires a closed surface. With the tip open, flow through the tip would 
cause erroneous results in this region and this was verified by the 
existence of significantly non-zero values of WNU (as computed from the 
velocity influence coefficients) at control points on the upper and 
lower panel columns adjacent to the tip. 

4.7 FLOW PAST ON ARROW WING-BODY CONFIGURATION 

An analysis was also performed on the arrow wing-body configuration 
shown in Fig. 2.5. This configuration has been the subject of an extensive 
series of wind tunnel tests ranging in Mach numbers from 0.40 to 2.50 
including supersonic Mach numbers as low as 1.05, reference 7 and 8. 

The configuration for which the comparisons are made featured a flat 
(no twist) wing with a rounded leading edge. The paneling representation 
is illustrated in Fig. 4.16 and 4.17. Details of the wing body inter- 
section are shown in Fig. 4.18 and 4.19. Composite source-doublet 
singularities are used to represent the configuration. The boundary 
condition imposed specify the source strength so as to cancel the normal 

component of the freestream mass flux (w, + Y ) • ¥ = 0 and the doublet 

' u 00 

distribution to cause the perturbation potential to vanish inside the 
configuration, 0^ - 0. Details of the boundary condition specifications 
were given in Section 2.5.4. 
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Figure 4.14 Test- Theory Comparison Wing Pressure 
Distribution Carlson Wing 2 
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A comparison of wing surface pressure distributions is shown in Fig. 4.20 
and 4.22, for a Mach number 1.70 and 2° angle of attack. Experimental 
data are compared to calculate results from the supersonic panel code 
and to 1 inearl ized supersonic results calculated by the FLEXSTAB programs, 
References 9 and 10. Note that the supersonic panel code results are 
taken to lie along the panel centerlines and do not necessarily correspond 
to the wing buttline at which the pressure measurements were made. The 
FLEXSTAB results, however, have been interpolated to the same butt! ires 
as the experimental data. The comparisons show general good agreement 
between the two theories and the experimental data. The lack of better 
test-theory agreement near the upper surface leading edge on the out- 
board portion of the wing is due to flow separation and the beginning 
of the formation of a leading edge vortex. Experimental isobar plots 
shown in Fig. 4.23 illustrate the presence of the leading edge vortex 
at this low angle of attack. This configuration is characterized by a 
very slender body and a thin highly swept wing which lends itself very 
well to the linearized approximation of the FLEXSTAB code. At first 
glance, both theories seem to yield very comparable results, but upon 
close examination it is clear the supersonic advanced panel pilot code 
is in better agreement with experiment near the leading 2 dge where the 
linearized approximation tends to breakdown. The predicted results 
from the supersonic panel code at the leading edge of the two inboard 
wing stations have been degraded by poor ecge matching at the wing body 
intersection. 

Fig. 4.24 compares surface pressure distributions along the body. Again, 
the location of the supersonic panel results do not exactly correspond 
to location of the experimental measurements. Very good agreement is 
seen between the two theories and the experimental data except for the 
supersonic panel results near the leading edge of the wing body inter- 
section. As previously stated, this descrepancy is due to a problem 
in the network edge matching in this area. 
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4.8 FLOW PAST A CONFIGURATION FEATURING SUPER-INCLINED PANELS 
An example of the use of superinclined panels is shown by the configuration 
in Fig. 4.25. Here an axi symmetric nacel le is modeled. A plane super- 
inclined network was placed inside a nacelle at axial position 2.25 as 
shown in the figure. This was done as a test case to cancel the distur- 
bance due to the upstream portion of the nacelle. In practice the super- 
inclined network could be placed at the face of the nacelle and would be 
necessary to cancel disturbance from the upstream portion of the remainder 
of the configuration. On the portion of the nacelle ahead of the super- 
inclined networks zero normal mass flux boundary conditions, w u • rf = 0 , 
were prescribed on the upper (or outer) surface of the source network. 
Potential boundary conditions (Morino-type) were prescribed on the 
remaining surfaces. Fig. 4.25 shows the exterior surface pressure distri- 
bution on the nacelle as a function of axial position for two different 
panel ings. The 48 panel network is seen to be too coarse to absorb 
adequately the interior reflections. Doubling the panels by increasing 
the number of radial divisions to form 96 panels smooths out the pressure 
distribution. The pressure is seen to be in good agreement with Li ghth ill's 
theory. 


i 


The results shown in this section have shown very reasonable comparisons 
between advanced panel pilot code and other theories and experimental data. 
The limited number of cases studied could in no way cover the spectrum of 
problems one may wish to analyze in the future. Nor, have any configura- 
tions of complex geometry which would better illustrate the advantages of 
the advanced panel pilot code over existing linearized theories been 
attempted. Never- the-1 ess, the results seen to date have been very 
encouraging as to the success of this new aerodynamic analysis method. 
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Appendix A - Singularity Splines 


SOURCE DISTRIBUTION 


On each panel, a plane may be defined which passes through the panel 
center and the panel edge midpoints. The source strength on a panel is 
taken to be linear in the coordinates of this plane: 

a = a o + a x X + 0 y Y 

The three parameters a 0 » a x > and o y which are needed to specify 
the source strength on each panel are not, per se, basic unknowns of the 
problem. Rather, they are calculated as needed in terms of the source 
strengths at the centers of the panels, which in turn are determined so 
as to meet the specified boundary conditions. 


Specifically, the parameters of Eq. (1A) are determined by a weighted 
least-squares fit to the source strength at the centers* of the eight 
surrounding panels. That is, if the source strength at the center of the 
ith of the 9 panels involved is called a -j , and the coordinates 
of the ith panel center are x-, y^ , the quantity 


Z »,(«, - X, - cr y Y 1 .) ; 


is made stationery with respect to variations m a cr , and cr 

a u * y 

The weight w = 10 H if i is the label of the center of the panel under 
study, and w = 1 otherwise. 


If the panel is at the edge of the network, there are fewer surrounding 
panels, but always enough to make this system determinate. Even when the 


*By "center" of a panel we mean just the average position of 
the four corner points. If two points are coincident, no special 
action is taken in locating the panel center, which thus weights 
the collapsed point twice as heavily as the other two. 


panel Is at a corner of the network, there are three neighboring panels, 
so the only effect is to reduce the number of terms in the sum which 
is minimized. See Fig. A.l. 

Two points are worth noting about this representation of the source 
strength. First, there is no reason for the source strength to be contin- 
uous from one panel to the next, since there is no special consideration 
given to the source strength at interpanel boundaries. Of course, in 
the limit as the grid is refined, we expcet the discontinuities in the 
source strength to vanish, so that these discontinuities furnish a 
measure of whether the mesh is sufficiently fine. Secondly, since the 
source distribution within the panels is determined in terms of the source 
strength at the centers of nearby panels, the number of parameters 
governing the source strength is exactly equal to the number of panels. 

DOUBLET STRENGTH 

The jump in normal mass flux (conormal velocity) across a surface singul- 
arity distribution is directly equal to the source strength. However, 
the jump in tangential velocity across a singularity distribution is 
equal to the gradient of the doublet strength. Hence it is necessary 
to use a higher-order representation of the doublet strength than of 
the source strength to assure the same level of representation of their 
velocity fields. Thus the doublet strength is represented by a quadratic: 

P - P 0 + ^ + ^y ^ + y xx 1 + ^xy ^ + *Vy 7 

where (x,y) are now the coordinates of points in a subpanel, a different 
quadratic being used for each subpanel. 
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O FIGURE A. 1 PANEL USED IN DETERMINING SOURCE STRENGTH 
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The determination of the doublet strength is complicated both by the necessity 
to use six parameters to specify the quadratic distribution and by the 
nesessity to avoid discontinuities in the doublet strength at either panel 
or subpanel boundaries. The determination is accomplished in stages, as 
follows: 

1. For each panel, nine "doublet-strength parameters" are defined, which 
may be associated with the doublet strength at the nine vertices of the eight 
triangular subpanels. As will be seen below in Fig. A. 2, 



FIGURE A. 2 SUBPANELS OF A DOUBLET PANEL 

these definitions are such that the doublet-strength parameters associated 
with points on the boundary between two adjacent panels are exactly the same, 
regardless of which panel is under study. This helps to achieve the desired 
continuity of doublet-strength across panel boundaries. 

2. Six of the ten straight lines which bound subpanels contain three 
subpanel vertices. In Fig. A. 2, these are the lines ABC, DEF, GHI, ADG, 

BEH, and CFI. The doublet-strength parameters associated with the vertices 
may therefore be fitted with quadratics along each of these six lines. 
Quadratics may also be used to describe the doublet-strength variation along 
the other four lines (BD, DH, HF, FB in Fig. A. 2), as follows. At point D, 
for example, the derivatives of the doublet-strength in the directions of 
the straight lines ADG and DEF may be computed by differentiating the quad- 
ratics already fitted to the doublet-strength parameters along those Tines. 
These data determine the derivatives of the doublet-strength at D along the 
lines DH and BD. Similar operations supply directional derivatives of the 
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doublet-strength at the other three midpoints of panel edges. A quadratic 
is then fitted to the doublet-strength variation along DH, using the values 
of the doubl et-strength parameters at D and H and, in a least-square fashion, 
the derivatives of the doublet-strength at those points. 

3. We now have quadratic representations for the doublet-strength 
along each subpanel boundary in terms of the doublet-strength at the nine 
vertices of the subpanels. These quadratics may be used to estimate 



FIGURE A. 3 SUBPANEL CONTROL POINTS 

the doublet-strength at the midpoints of the subpanel edges (X's in Fig. 

A. 3) in terms of the same nine parameters. Then, for each subpanel, we have 
six bits of information on the doublet-strength; its valuesat the three 
vertices of the subpanel and at the midpoints of its three sides. These 
six bits of data are just enough to determine the six coefficients of a 
quadratic representation of the doublet-strength within the subpanel. 

4. Thus, within each panel, we have a piecewise-quadratic representation 
of the doublet-strength in terms of nine doublet-strength parameters. The 
representation can be seen to be continuous across subpanel boundaries, as 
follows. The quadratics in the neighboring subpanels BCF and BEF (see Fig. 
A. 4) are each determined, in part, by values of the doublet 


E F 



FIGURE A.4 SUB PANEL CONTROL POINTS' 
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strength at B, F, and the midpoint of BF, Along FG, these quadratics be- 
come quadratic in the distance along BF, Since a quadratic in one variable 
is determined uniquely by three values, these quadratics must agree at every 
point along BF, not just at the three interpolation points. Similarly, since 
the singularity parameters associated with points along panel edges are the 
same for the panels on either side of the edge, the piecewise quadratic repre- 
sentation is continuous between panels. However, its normal derivative may 
be discontinuous both at panel and subpanel boundaries. 



5. It remains to define the nine doublet-strength parameters. As in the 
case of the source distribution, it is convenient to regard as basic unknowns 
the doublet-strength at the panel centers. This helps to assure a close 
corresponderce between the number of unknowns and the number of panels. For 
panels well within the network interior, the doublet-strength parameter at 
a corner of a panel is determined by least-squares fitting of a quadratic to 
the doublet-strength at the centers of 12 nearby panels, as shown in Fig. A. 5a. 
Again, the squared errors at the four panel midpoints closest to the corner 
point in question are weighted 10® times as much as the squared errors at the 
outer six panel midpoints. Similarly, the doublet-strength parameters at the 
midpoint of a panel edge is found by least-squares fitting a quadratic to the 
doublet strength at the center of eight nearby panels, the squared errors at 
the centers of the two panels bordering the edge being weighted 10® more 
heavily than the other errors. Finally, the doublet-strength parameter at 
the center of a panel is simply set equal to the doublet-strength at that 
point. These arrangements are illustrated in Fig. A 5b and A. 5c. When a 
panel corner point or edge midpoint is near the network boundary, the fixing 
of the doublet-strength parameter at that point can be made similar to the 
procedure employed in the interior by introducing, as additional basic un- 
knowns, quantities related to the doublet-strength at the midpoints of the 
panel edges which border the network. Fig. A 6a through A. 6e illustrate 
the utility of this arrangement. We also add to our list of basic unknowns 
the doublet-strength at the corner of the network, which then becomes the 
doublet strength parameter for the corresponding corner of the corner panel; 
see Fig. A. 6f. 
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Then the doublet-strength parameters at the network corner and at the mid- 
point of a panel edge which lies on the network edge are determined directly 
in terms of the corresponding unknown; see Fig. A.6f and A.6g. Finally, the 
parameter at a panel corner point which lies on the network edge is determined 
by least-squares fitting a quadratic to the doublet-strength at the four neigh- 
boring point shown in Fig. A.6h. 

The upshot is that the number of unknowns which govern the doublet-strength 
equals the number of panels, plus the number of panel edges which border 
the network, plus the number of network corner points. 

Naturally, special attention must be given in cases for which the panels, as 
illustrated in Fig. A. 7, are triangular rather than quadrilateral, if only 
because it is impossible to fit a quadratic to three values at a common point! 
However, the same basic selection of unknowns serves to define enough doub- 
let-strength parameters to fix piecewise continuous quadratics in each of the 
six subpanels of any triangular panel. 



FIGURE A. 7 SUBPANELS OF TRIANGULAR PANEL 
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